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ABSTRACT 

We present one of the most precise measurement to date of the spatial clustering of X-ray selected AGNs 
using a sample derived from the Chandra X-ray Observatory survey in the Bootes field. The real-space two- 
point correlation function over a redshift interval from z = 0.17 to z ~ 3 is well described by the power 
law, £(r) = (r/ro)~ r , for comoving separations r < 20h~ l Mpc. We find y = 1.84 ± 0.12 and ro consis- 
tent with no redshift trend within the sample (varying between ro = 5.5 ± 0.6 h Mpc for (z) = 0.37 and 
ro = 6.9 + 1.0 h Mpc for (z) = 1.28). Further, we are able to measure the projections of the two-point 
correlation function both on the sky plane and in the line of sight. We use these measurements to show that the 
Chandra/Bootes AGNs are predominantly located at the centers of dark matter halos with the circular velocity 
Vmax > 320 km s _1 or M\%q > 4.1 x 10 12 /i~' M Q , and tend to avoid satellite galaxies in halos of this or higher 
mass. The halo occupation properties inferred from the clustering properties of Chandra/Bootes AGNs — the 
mass scale of the parent dark matter halos, the lack of signific ant redshift e volution of the clustering length, 
and the low satellite fraction — are broadly consistent with the Hopkins et al. (2006 ) scenario of quasar activity 
triggered by mergers of similarly-sized galaxies. 



1. INTRODUCTION 

Direct observations of host galaxies of high-redshift active 
galactic nuclei are hard or impossible with the current instru- 
mentation, except for highly obscured or low-luminosity ob- 
jects. Hence, studies of the AGN clustering properties are a 
unique source of information on the AGN hosts and their envi- 
ronment. At low redshifts, the supermassive black holes exist 
at the centers of most low-redshift galaxies, and there is a tight 
correlation between the SMBH mas s and the properties of the 
bulges of host galaxies (jM agorrian ~et al.||1998| Ferrarese & 
Ford 2005 Gebhardt et al.||2000| >. This suggests that most 



galaxies hosted an AGN at some point in their evolution, and 
that AGNs at each redshift are stochastic "markers" of a pop- 
ulation of galaxies in which the conditions are favorable for 
accretion of matter onto the central SMBH. Through match- 
ing the clustering properties of AGNs to those of dark mat- 
ter halos, or with those of different types of galaxies, or with 
those of AGNs of different types, we can determine the typ- 
ical mass scale of the AGN hosts, their morphological type, 
and determine whether different types of AGNs are hosted in 
the same type of objects. 
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Clustering of optical quasars indeed is very similar to that 
of galaxies. The two-point correlation function observed at 
z ~ 0.5 - 2.5 in the separation range 1 - 20h~ l Mpc is well 
described by a power law, £(r) = (r/>o)~ y with a slo pe of 
y — 1.9 and a correlation len gth of ro = 5.5 h~ l Mpc (Ross 
iraLl[2009l |Shen et al.|2009). Further , usi ng the SDSSand 



2QZ quasar sample, Croo m et aL| ( |2005[ > and |Ross et al.| ( |2009 1 
track the evolution of the optical quasar clustering amplitude 
over the redshift range z = 0.5 - 2.5 and find only a mild 
evolution. 

The measurements of the spatial clustering of X-ray AGNs 



start w ith the works of Mullis et al. ( 2004 1 and Yang et al. 
(2006). The most recent measurements of the spatial autocor- 



relation function of X-ray AG Ns can be found in |Gilli et aL 
( [2003] [2009] l and |Cappelluti et al.| | |20l0) . In several recent 
works, the X-ray AGN clustering has beed studied through 



cross-correlation with galaxy catalogs (Hickox et al. 2009 



Coil et al.||2009| |Krumpe et~aT"]|2010|i. Generally, the two 



point correlation function of the X-ray AGNs was found to be 
similar to that of the optical quasars; however, the previous X- 
ray studies could not constrain the evolution of the correlation 
function over a sufficiently wide redshift range. 

The main conclusion from previous clustering analyses is 
that AGNs are located in gala xy group-sized dark matter ha- 
los (M ~ 2 x 10 12 ft -1 M , see |Ross et al.|2009| , with the mass 
scale fairl y independent of the object redshift or observed lu- 
minosity. Hick ox et al.| f2009 ) also find some differences in 
the clustering properties and color of host galaxies for the X- 
ray, radio, and infrared-selected AGNs at z = 0.5. The differ- 
ence in the AGN clustering properties and colors of their host 
galaxies lead 



techniques se 



Hickox et al. 



to conclude that these different 
lect distinct source populations and not simply 
different stages of rapidly changing AGN properties. 

In addition to determining the mass scale of the AGN host 
dark matter halos, it is interesting to establish where within 
the halos the active galaxies are located. This question can 
be addressed by direct observations of individual objects only 
in rich galaxy clusters. Indeed, some studies ( |Cappi et al. 
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Fig. 1. — The two-dimensional correlation function f(r p ,;r) for objects lo- 
cated within dark matter halos with the maximum circular velocity V max > 
250 km s at z = 1, plotted as a function of transverse (r p ) and line-of-sight 
(n) pair separation. 70% of objects are forced to be at the centers of such 
halos, while 30% are put in one of the satellite subhalos. In the left panel, the 
peculiar velocities were set to 0. 



[200T] |Molnar et aLl|2002l |Martini et aLl |2006) indicate an 
exce ss of X-ray AG Ns in the cluster outskirts (see, how- 
ever, [Koulouridis & Plionis 2010). However, the majority of 
quasars and AGNs are located within galaxy group-sized ob- 
jects, where it is hard or impossible to independently localize 
the centroid of the system, especially at high redshifts. For 
such systems, it is possible to determine the fraction of objects 
in satellite galaxies through analysis of the two-point correla- 
tion function at small separatio ns. An example of such an ap- 
proach can be found in Padmanabhan et al. (2009). Based on 
cross-correlation of the SDSS quasar and Large Red Galaxies 
samples, these authors argue that a high fraction, > 25%, of 
the optical quasars must be located in non-central galaxies. 

One of the main goals of the present work is to constrain 
the location of X-ray emitting AGNs within their host dark 
matter halos by a more direct method. The effect we are using 
is the strong dependence of a galaxy's peculiar velocity on its 
location within the host halo. The central galaxies are pre- 
dicted (and observed, at least in massive clusters, see Oegerle 
|& Hill|2001| l to be nearly at rest with respect to the host halo, 
and their random motions correspond to motions of the ha- 
los as a whole. The satellite galaxies move at approximately 
the virial velocities within their parent halo. As a result, the 
satellite galaxies have much faster peculiar motions and form 
"finger of God" structures in the radial velocity space. The ef- 
fect of high peculiar velocities of the satellite galaxies can be 
detected through comparison of the objects' clustering prop- 
erties as a function of projected separation, r p , and the line- 
of-sight separation, n. Peculiar motions do not affect r p but 
can strongly distort n because it is derived from the measured 
redshift. 

Using peculiar motions in the previous studies of the AGN 
clustering was hard because of the small number of objects 
and insufficient accuracy of the redshift measurement^ 2 ] For- 
tunately, we now can use an excellent sample for sucn stud- 
ies. The Chandra survey of 9.3 deg 2 in the Bootes region 
( Murr ayet al.|20 05 ) provides a uniformly selected sample of 
> 3000 X-ray selected AGNs, for a uniform subsample of 
which (~ 1900 sources) high-quality spectra were measured 
with MMT/Hectospec (Kochanek et al., in preparation). 

Full information on the object clustering and peculiar veloc- 

12 For example, the redshift uncertainties corresponded to peculiar veloci- 
ties of a 420 - 500 bus' 1 inthe2 QZ sample jCroom et al.|2005) . With such 
large uncertainties, Croo m et al.| could not extract useful information from 
the velocity-space distrortions although they were included in the modeling 
of the spatial correlation function. 




Fig. 2. — Chandra survey in the Bootes field. Grey scale map shows the 
combined Chandra sensitivity map (exposure maps for 1 26 individual point- 
ings, each multiplied by the poi nt so urce detection sensitivity as a function 
of off-axis distance derived in § |2.1( . Crosses mark the location of X-ray 
sources used in this work. The polygon shows conservative boundaries of 
the spectroscopic survey. We consider only 1282 sources within the polygon 
boundaries in the clustering analysis. 



ities is contained in the two-dimensional correlation function, 
g(r p , n). An example derived from the numerical simu lations 
we use in t his Paper is shown in Fig. [Tl See, e.g., [Zehavi 
et al. (2002} for an example of the detailed modeling of the 



f(r p ,7r) function measured for the SDSS galaxies. Unfortu- 
nately, such a detailed modeling is impossible for the present 
high-z AGN samples due to limited statistics. However, we 
show that useful information can still be obtained using pro- 
jections of ^(r p ,7r) on the radial velocity direction and on the 
sky plan e, %„{n) and £ p (r p ), respectively (formally defined 
in § |3.1| below). Using numerical cosmological simulations, 
we show that if a substantial fraction of objects lie in satel- 
lite galaxies, %„{n) is expected to be significantly in excess 
of £p(r v ) in the range of comoving distances 1 - 10 1 Mpc 
(§ |4.3| >. The obse rved correlation functions do not show such 
an excess (§ 5.1 1, which we exploit to put an upper limit on 
the fraction of X-ray AGNs in the satellite galaxies. Finally, 
we explore the redshift evolution of the clustering length and 
compute the typical mass scale of the AGN host dark matter 
halos and put constraints on the AGN duty cycle. 

All cosmology-dependent quantities are computed assum- 
ing a spatially fiat model with parameters Q M = 0.268 and 
Qa = 0.732 (best-fit ACDM paramaters obtained from the 
combination of CM B, superno vae, BAO, and galaxy cluster 
data, see Vikhlinin et al. 2009). All distances are comoving 
and given with explicit //-scaling, where the Hubble constant 
is Hq = 100 h km s Mpc 1 . The parameter uncertainties 
are quoted at a confidence level of 68%. 

2. AGN SAMPLE 

We use a sample of high-z AGNs derived from the Chan- 
dra X-ray survey in the 9.3 deg 2 Bootes field of the NOAO 
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Fig. 3. — The redshift-luminosity diagram for Chandra sources in the Bootes 
field used in our clustering analysis. 



Deep Wide-Field Survey (Murray et al. 2005). The region was 
uniformly covered with a grid of overlapping 5 ksec ACIS-I 



pointings providing a sensitivity of 4.7 x 10 15 erg s 1 cm 



-2 



in the 0.5-2 keV energy band (Ken ter et al.|2005] ). Extensive 
optical data exist for this field and 98% of X-ray sources have 
optical or infrared counterparts ( Brand et al.|2006| . The red- 
shifts for X-ray sources with optical counterparts brighter than 
I = 21.5 were uniformly obtained with the MMT/Hectospec 
instrument in the AGN and Galaxy Evolution Survey (AGES, 
C. S. Kochanek et al., in preparation). Further details on the 



X-ray and optical observations can be found in Hickox et al. 
( |2009l ). 

For clustering analysis, we use 1282 X-ray selected AGNs 
with spectroscopic redshifts. Redshift measurements effec- 
tively introduce an additional selection criterion for the sam- 
ple, based on the optical magnitude of the counterparts, / < 
21.5. This selection removes ~ 35% of the X-ray sources 
and likely introduces a hi gh-z cutoff in the redshift distribu- 
tion of our sources (§ |2.2| i. Fortuna tely, the redshift measu re- 
ments were done with Hectospec (Fabrican t et al.| |2005) in 
several passes, which excludes the so called "fiber collisions" 
problem. The locations of X-ray sources with spectroscopic 
redshifts are marked in Fig. [2] and their luminosity-redshift 
diagram is shown in Fig. [3] 

To reconstruct the correlation function of sources, we need 
to simulate catalogs of randomly distributed sources whose 
distribution follows the sensitivity variations and the redshift 
distribution of our catalog. Below, we describe how these 
functions were derived from our Chandra sample. 

2.1. Spatial variations of sensitivity 

Chandra sensitivity for the detection of point X-ray sources 
is not uniform across the field of view. Sensitivity variations 
(~ ±25%, see below) are imprinted on the distribution of de- 
tected sources. The typical spatial scale of the sensitivity vari- 
ations is several arcmin, which corresponds to the comoving 
distance of 1 - 2 hr x Mpc at z — 1. This is comparable to the 
distance scale where we measure clustering, therefore these 
variations must be taken into account. 

The sensitivity variations are caused mainly by two effects 
— 1) the vignetting of the Chandra X-ray telescopes, and 
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Fig. 4. — The observed surface density of detected X-ray sources as a func- 
tion of the distance from the telescope optical axis. The sharp decrease at 
r c < 1' is attributed to the gap between the Chandra ACIS-I CCDs. The 
gradual decline at larger radii reflects the sensitivity variation caused by mir- 
ror vignetting and PSF degradation. The solid line is the best-fit second-order 
polynomial (eq. [T]). 

2) the degradation of the Point Spread Function (PSF) away 
from the optical axis. The mirror vignetting is well-calibrated 
and its effect can be computed for a given source population. 
However, the effects of PSF degradation on the source detec- 
tion efficiency in the low-photon regime are very complex. 
Therefore, it is best to measure the combined effect of the 
sensitivity variations empirically. 

Both the mirror vignetting and PSF degradation are approx- 
imately azimuthally symmetric. Therefore, we can assume 
that the detection sensitivity is a function of the source offaxis 
angle, r c and measure it from the radial profile of the surface 
density distribution of detected X-ray sources, averaged over 
all the 126 fields. The results are shown in Fig. HI For a 
uniform sensitivity, we would expect a constant surface den- 
sity, while in reality we find that the average surface density is 
~ 7 arcmin -2 near the optical axis, and falls to ~ 5 arcmin 2 
at off-axis distances of 10 arcmin. The derived radial profile 
of the source surface density can be modeled with a second- 
order polynomial, 



dN 
ds 



= 6.9(l-(r f /17.4') 2 ). 



(1) 



In addition to this gradual variation with radius, there are 
sharp features in the spatial sensitivity pattern related to the 
gaps between the ACIS-I CCDs. In particular, these gaps are 
responsible for the drop in the number of detected sources at 
r c < 1' in Fig. |4j These sensitivity variations can be ade- 
quately taken into account using the standard Chandra expo- 
sure maps. In doing so, we assume that the expected surface 
density of X-ray sources at a given location is proportional 
to the effective exposure at this location. This assumption is 
justified because the source detection in the Chandra/Bootes 
survey is well in the photon-limited regime (the background is 
unimportant) and the observed cumulative source counts are 
very close to N(> S) oc S" 1 ar ound our flux limits ( Hasinger 
|et al.|1993HVikhlinin et al.jl995] l. 

Our final sensitivity map (Fig. [2|l consists of the merged 
set of Chandra exposure maps computed for the individual 
126 pointings, each multiplied by the radial sensitivity pattern 
given by eq. (QJ. This sensitivity map is taken into account in 
the derivation of the correlation function through the genera- 
tion of the appropriate catalogs of random sources. The total 
area within the conservative boundaries of the spectroscopi- 
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cally surveyed region of the Bootes field is 7.30 deg 2 . The 
total effective area, taking into account the gaps between the 
Chandra CCDs and degradation of the detection efficiency, is 
5.90 deg 2 . 

2.2. Redshift distribution 

The model for the source redshift distribution, dN/dz, 
should reflect both the intrinsic variations of the comoving 
number density with redshift and all selection effects of the 
catalog. A commonly used approach is to model the observed 
dN/dz distribution with a high-order polynomial (see, e.g., 



Croom et al.||2005| |Ross et aL]|2009| for recent examples) 
This approach works well for catalogs with a large number 
of sources. However, for smaller catalogs, like ours, there is 
a danger that a high-order polynomial fit will follow statis- 
tical fluctuations in the observed dN/dz, while a low -order 
polynomial would be unable to adequately model the strong 
gradients at low z. Therefore, we fit the redshift distribution 
of AGNs in the Bootes field with a parametric model based 
on several physical assumptions. 

The first component of the model represents the cosmologi- 
cal comoving volume per unit redshift, dN\/dz °c dV/dz- The 
second component is a power law function of the minimum 
luminosity which corresponds to the Chandra flux limit at 
redshift z, dNi/dz °c L a J^ n cc d", where di is the luminosity 
distance. This component represents the effect of the low-L x 
cutoff of the intrinsic luminosity function introduced by the 
selection which is primarily based on Chandra detections. It 
also can describe the evolution of the luminosity function at 
high z. The third component is a high-z cutoff modeled by 

a broad Gaussian, dN^/dz °c exp [~d 2 L /C 2 ^. This component 
can represent the high-L A cutoff or steepening of the intrin- 
sic AGN luminosity function, and also can describe various 
observational limits implicitly built into our catalog (e.g., a 
lower efficiency of optical identifications and redshift mea- 
surements for the highest-z X-ray sources). This simple ana- 
lytic model, 



dN 

~dz 



dV 

dz 



= consl x — d a L exp(-(<f L /C) 2 ), 



(2) 



which has only two free parameters provides a strikingly good 
fit to the observed redshift distribution of the Bootes X-ray 
selected AGNs (Fig.lSj. The best-fit values are a — -1.07 and 
C = 1.50x 10 3 h~ l Mpc. We use the arguments outlined above 
only as a motivation for a good analytical description of the 
dN/dz distribution for our sources. The functional form and 
derived parameters are not meant to represent the true X-ray 
luminosity function or its evolution. 

Figureplalso demonstrates the general characteristics of our 
sample. The peak in the observed dN/dz distribution is near 
z m 0.6. The median redshift of the sample is z me d = 1.04. 
The tail in the redshift distribution extends to z ~ 4.5 but the 
fraction of AGNs with z > 3 is very small. Overall, the clus- 
tering properties of sources in our sample are most sensitive 
to the distribution of the X-ray AGN population near z ~ 1. 

3 . TWO-POINT CORRELATION FUNCTION OF Chandra/Bootes AGNS 

3.1. Definitions 

Because the volume is never sampled completely in astro- 
nomical surveys, the derivation of the two-point correlation 
function from the data uses mock catalogs of intrinsically ran- 
domly distributed objects, which faithfully reproduces all ob- 
servational distortions introduced by the survey. Examples of 



such distortions are boundaries of the survey region, gaps in 
the data or spatial variations of the sensitivity, variations of 
the selection efficiency with redshift, etc. Given the catalog 
of observed sources and the mock random catalog, the two- 
point correlation function can be estimated (Landy & Szalay 
^993) as 



N r (N t -I) DP N T -l PR 
N d (N d -1)~RR ~ N d RR 



+ 1, 



(3) 



where DD is the number of source pairs in the data for the 
given distance interval, RR is the corresponding number of 
pairs in the random catalog, DR is the number of pairs be- 
tween the data and random catalog, and Na and N r are the 
numbers of objects in the data and random catalogs, respec- 
tively. Statistical uncertainties for £ can be estimated as 



(4) 



(Peebles 1973[ ); this equation includes both the Poissonian 
shot noise and intrinsic variance terms. To verify the accu- 
racy of the error by eq.RJ we used the sample varience of the 
correlation functions measured in t he mock catalogs derived 
from the Millenium simulation (Kitzbi chler & White||2007] l 
for the survey geometry and object properties similar to those 
of the Chandra/Bootes survey. This analysis showed that eq.|4] 
is accurate at tsmall scales but may underestimate the uncer- 
tainties at large scales. The correction factor can be described 
by a smooth function which is 3%, 23%, and 42% at separa- 
tions of 1, 6, and 15 ft -1 Mpc, respectively. This correction is 
applied to the statistical uncertainties estimated by eq. [4] 

The correlation function in real space is expectedto be 
isotropic, so £ is a function of the 3D separation only. When 
the object redshifts are used to derive the distances, the cor- 
relation function is distorted in the line-of-sight direction be- 
cause of large-scale flows (the Kaiser |1987| effect) and "fin- 
gers of God" arising within the virialized dark matter halos. 
The correlation function should then be measured as a func- 
tion of the projected separation, r p , and the line-of-sight sepa- 
ration, 7r. Equations (3) and HI still can be used, but the pairs 
must be counted for each combination (r p , ri). 

Given the angular separation between two objects, 8, and 
redshifts, zi and zi, the comoving separations r p and n can 
be computed as follows. First, one computes the radial co- 
moving distances, D c \ and D c ,i , corresponding to t he object 



redshifts (see, e.g., Hogg 1999 1. Then, following Davis & 
|Peebles| ( fT983] l we have 



\D c ,i - D c 



r p = [2D Ci iD c , 2 (l-cos0)] 



1/2 



One can also define a formal 3D separation, 

a/2 



(5) 



(6) 



(7) 



but it should be kept in mind that s is not equivalent to the true 
3D separation, r, because of the redshift space distortions. 

As only the line-of-sight separations, n, are affected by the 
object peculiar velocities, it is useful to consider the correla- 
tion function projected on the sky plane, 



w p {r t 



(8) 
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Fig. 5. — The observed distribution of Chandra/Bootes AGNs as a function of redshift (left) and comoving distance (right). The best fit by eq. J2] is shown by 
solid lines. 

this analysis, we trancated the integration at comoving sep- 
arations ±30 /T 1 Mpc, and further reduced the sensitivity of 
our results to this choice by applying small corrections to the 
integrated functions. Assuming that the true correlation func- 
tion follows a power law with index y at large separations, we 
can calculate the effect of this truncation as 



because it is not mo dified by the redshift-space distortions 
( |Davis&Peebles|1983) : 

£ £ (tme) ( yjrl+K 2 ) dn = £ e hs) (r p , n) dn. (9) 

Because of the insensitivity of w p (r p ) to the redshift-space 
distortions, most of the studies which involve detailed mod- 
eling of the shape of the galaxy two-point correlatio n func- 
tion are bas ed on fitting the w p (/p) mea surements (e.g., |Zehav7 
et al.|2005 Padmanabha n et al.|2009 1. Using w p (r p ) instead 
of f (r) is particularly straightforward in those cases when £(r) 
can be sufficiently accurately approximated by a power law, 
£(r) = (r/ro)~ r . It follows from eq. |8) that in this ca se the 
relation between £ and w p (r p ) is simply ( |Peebles| 1980) 



where 



w p (r p ) = A(y) r p (r p /r ) 



A(y)=r(l/2)n(y-l)/2]/T(y/2). 



(10) 



(11) 



Therefore, the correlation length, ro, and the slope of the true 
3D correlation function £ can be obtained immediately from 
the power-law fit to w p (r p ). 

Below, we also use a projection of the measured 3D corre- 
lation function on the line-of-sight direction, 



w„(n) = I 

J -ex 



? ohs \r p ,n)dr p . 



(12) 



In the absence of peculiar motions, w„(n) should be equivalent 
to w p (r p ). In particular, for a power-law £(r), eq. |T0) will be 
valid also for w„(ri). 

Large separations do not contribute significantly to the in- 
tegrals in equations [8] and [9] but add noise, therefore in prac- 
tice the integration is truncated at some finite separation. A 
good choice for studies of AGN sample s is to truncate th e 



integration at as +40/; Mpc (see, e.g., |Gilli et al. 



2009) . 



Following this and other works we investigated how the de- 
rived correlation length depe nds on the cho ice of the cutoff 
distance (see, e. g., Fig. 9 in |Gilli et al.|20 05 or Fig. 3 in Alle- 
vato et al.|20lT i. For our sample, the convergence is reached 
at a cutoff radius of 20 Mpc (for cutoffs between 20 and 
60h~ l Mpc, the derived ro varies by less than 4%). Based on 



C y {r) = 



-7/2 



dl 



gp + p-yf* d i 



(13) 



The correction coefficient, C y (r), is close to 1 for r «: 
30/i _I Mpc and gradually decreases at larger distances. For 
y -2 (close to our best-fit value), C y = 0.96 at r — 2 h Mpc, 
0.87 at r — blT 1 Mpc (w the correlation length), and 0.63 at 
r = 20/i~'Mpc (« the largest separation at which the pro- 
jected correlation is still marginally detectable). 

Using equations 10 and 13 we can define the quantities 



1 



w p (r p ) 



A(y)CJr p ) r 



1 



A(y)C y (n) n 



(14) 



(15) 



where y is determined from the power-law fit to the measured 
w p (r p ) (eq. p0|). Thus defined, % p (r p ) should be a close ap- 
proximation to the true 3D correlation function ^(r). In the ab- 
sence of redshift-space distortions, ^„(n) also should approx- 
imate £(r). However, the peculiar motions (and in particular, 
the "finger of God" effect) suppress ^„{n) on the smallest sep- 
arations and enhance it relative to %(r) on intermediate scales 
(see below). Therefore, % p (r v ) can be used to determine the 
correlation length of our objects and hence the mass scale of 
their host dark matter halos. The ratio ^„{7i)l^ p {r v ) reflects 
the amplitude of the peculiar motions and hence can be used 
to constrain the fraction of objects located in the satellite dark 
matter subhalos. 

3.2. Results 

The correlation functions £ p (r p ) and ^„(tt) derived for our 
complete sample are shown in Fig. [6^ (filled and open cir- 
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Fig. 6. — (a) The two-point correl ation functions of Bootes/Chandra AGNs estimated from projections on the sky plane (fp(/p), filled circles) and on the line 
of sight (fjr(7r), open circles). See § |3.2| for the description of our procedure to derive fp(r p ) and t; n (n). Dashed line shows a power law fit to the f p (r p ) data, (b) 
The ratio of two projected functions. At small separation (n < 1 h Mpc, indicated by the vertical dashed line), the ratio can be affected by uncertainties of the 
AGES redshift measurements. At larger separations, we observe no statistically significant difference between t; p (r v ) and 



cles, respectively). Technically, we derive the functions as fol- 
lows. First, we estimate the two-dimensional correlation func- 
tion, £;{r p ,ji), on a grid of separations with equal logarithmic 
width for each cell, Alogc/ = 0.2. We use the Landy & Sza- 
lay estimator (eq. BJ. The random catalog is generated us- 
ing the spatial sensitivity ma p and the mo del for the redshift 
distribution in our sample (§ 2.1 and 2.2 1. To minimize the 
additional noise, the number of objects in the random cata- 
log is a factor of 100 larger than that for the AGN sample. 
The derived £(r p ,n) is then integrated in the radial and sky 
plane directions to obtain w p (r p ) and w n {ji) (eq. [9] and [l2j i 



As discussed in § 
tions 



3.1 



the integration is limited to separa- 
30 hT l Mpc to minimize noise. We then fit a power 
law (eq. 10 1 to w p (r p ) to measure the slope of the correlation 
function, y = 1.84 + 0.12. The correlation length measured 
for the full sample is r = 6.41 ± 0.44 h' 1 Mpc. We then 
compute the renormalization factor A(y) (eq. 11 1 and correct 



for the truncation of the integration in W p {rp) and w„{n) at 
1 Mpc (eq. 



13 I. With the renormalization factors deter- 



30 hr 

mined, we convert w p (r p ) and w„(n) into % p (r p ) and ^(n) us- 
ingeq. (T4][T5). 

As discussed above, £ p (r p ) should be close to the true 3D 
correlation function %(r), while %„{n) should be distorted by 
peculiar motions of the AGN host galaxies. Indeed, we ob- 
serve a suppression in ^{n) relative to % p (r p ) at separations 
n < Ih Mpc. Unfortunately, the corresponding radial ve- 
locity difference, Az m 0.0003, is uncomfortably close to the 
uncertainties in the AGES redshift measurements (Kochanek 
et al., in preparation). Therefore, we ignore the £„(tt) data 
at n < 1 /i Mpc. At intermediate separations (r » [1 — 
10] hT l Mpc), £;„{ji) is expected to be enhanced by the pecu- 
liar motions (see below). No such enhancement is present in 
the data. In fact, the ratio £n(?f)/£p( r p) I s fully consistent with 
1 at separations > 1 h~ l Mpc within the measurement uncer- 
tainties (Fig. j^i). Below, we show that this can be used to 

13 The uncertainties of the f;r(X)/£p('"p) ratio are computed by propagation 
of uncertainties in (;„(n) and t; p (r v ). In principle, a region in the 2-dimentional 
i;{r p ,n) enters into both projections. However, the analysis shows that the 



put an upper limit on the fraction of AGNs that can reside in 
satellite galaxies orbiting within massive dark matter halos. 

3.2.1. Comparison with previous measurements 

Our results represent the most accurate measurement of the 
spatial clustering of X-ray selected AGNs to date, so a com- 
parison with earlier observations is useful. 

The first detection of angular cluster i ng of X-ray sources 
was reported by Vikhlinin & Forman ( \ 995) based on the 
analysis of the ROSAT PSPC data. Using the Limber equa- 
tion reconstruction (Peebles 1980), these authors estimated a 
correlation length of 19 ± 5 h~ L Mpc (comoving) at z = 1.5 
using raw measurements, and 7.5 ± 2hT l Mpc after correct- 
ing for the "amplification" bias caused by the poor angular 
resolution of the ROSAT PSPC. The latter value is in good 
agreement with our results. 

Our results also are in good agreement with direct mea- 
surem ents of the spatial c lustering in the ROSAT NEP sur- 
vey ( jCappellutTe t al. 200 7|), Chandr a surveys in Deep Fields 
North and Sout h (|Gilli et al.||2005) , and the XMM-Newton 
COSMOS field ( |Gilli et al |2009[ >. The correlation length we 
measure in the low redshift bins (see 5 5.2 below) is in good 
agreement with rp = 5.5 - 6/T 1 Mpc found in two studies 
based on cross-correlations of X-ra y AGNs with galaxy cata- 
logs |Hickox et al.| ( |2009) ; |Coil et al^2009[ l, although is higher 
than the rp = 4.3 + 0.35 h Mpc reported in a similar work by 



Krumpe et al. 

Correlation lengths for the X-ray AGNs have been esti- 
mated also from the Limber inversions of the angular clus- 
tering measured for the XMM-Newton and Chandra sources. 
A wide range of rp values can be found in the literature (e.g., 
|Plionis et al |2008| and references therein); our measurements 
are inconsistent with rp > 10 hr x Mpc reported in some of 
these analyses. 

4. AGN CLUSTERING MODEL 



data in this region contribute negligibly to the "signal" and "noise" in the 
ratio, so our simple uncertainty calculation is adequate. 
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In the past ten years, a very successful framework for mod- 
eling the nonlinear clustering properties of galaxies has been 
developed ( the so called Halo O ccupation Distribution (HOD) 
model, see |Seljak|[2000l IMa & Fry|[2000l IPeacock & Smith l 
200 01 |White et al.|200ll|5cb ccimarro et al.|2001[|Berlind &| 
Weinberg 2002; Berlind et al. 2003 among others). The ap- 
proach is based on the idea that the distribution of dark mat- 
ter can be fully described through the mass function, linear 
bias, and density profiles of dark matter halos. These elements 
are well-calibrated using A^-body simulations. The two addi- 
tional ingredients of the model, which are less well known, 
are the probability distribution for a halo of mass M to con- 
tain N galaxies, and the distribution of galaxies within the 
halos. These functions can be parameterized by the functions 
suggested by the results of high-resolution numerical simula- 
tions (e.g., Kravtsov et al. 2004), and some parameters of the 
model can be in fact determined by fitting the observe d cor- 
relation func t ions. In particular, Kravtsov et al. (2004| and 
Zheng et al. ( |2005 1 show that the elements of HOD can be 
effectively decomposed into two components, separately de- 
scribing the properties of central and satellite galaxies within 
the dark matter halos. 

The HOD model is now very well developed for fitting the 
shape of the two-point correlation function. This technique 
has been applied to modeling the projected two-point correla- 
tion functions, w p (r p ), for objects ra nging from Lyman-break 
galaxies at z = 3 - 5 ( jConroy eT al. 2006) to relatively low- 
z quasars from SDSS (Padmanabhan et al. 2009). Recently, 
the HOD approach has been developed also for modeling the 
redshift-space distortions in th e galaxy correlation functions 
( Tinker et al.|2006 Tinker|20"07 1 — just the type of informa- 
tion we are aiming to use in this work to constrain the loca- 
tions of Bootes/Chandra AGNs within the host dark matter 
halos. 

In principle, the HOD models for galaxy clustering are an- 
alytic, and thus are convenient for those applications in which 
the cosmological parameters are varied. However, some of the 
most essential parameters of the HOD models are calibrated 
using numerical simulations. If one is interested in varying 
the parameters of galaxy distribution at a fixed redshift in a 
fixed cosmological model, it is more accurate — and easier — 
to obtain the model correlation functions directly from numer- 
ical simulations rather than to rely on analytic approximations 
derived from analyzing the simulations. This is the approach 
we take here. 

4.1. Numerical Simulations 

The set of simulations we use in this work is described 
Tasitsiomi et al.| p004[ ) and |Conroy et al.| ( |2006| ). These 



in 

are high-resolution dissipationless simulations run in a flat 
ACDM cosmology with parameters close to the present-day 
"concordance" values, Qm = 0.3, h = 0.7, and <x 8 = 0.9. 
The simulations follow the evolution of dark matter in a 
120 h~ l Mpc box. The box contained 512 3 dark matter par- 
ticles with mass m p = 1.07 x 10 9 hr l M Q ; the peak resolution 
reached 1.8 hr l kpc. 

The locations and velocities of the dark matter particles 
in the simulations were recorded at z = 0.09, 0.5, 1.0, 2.0, 
3.1, 4.0, and 5.0, matching well the redshift distribution in 
the Bootes/Chandra AGN sample (see Fig|5]l. The simula- 
tion outputs were then analyzed to identify dar k matter ha- 
los an d subhalos using a modification of the Klypin et al. 
( |1999[ ) bound density maxima halo-finding algorithm (see 



Kravtsov et al.|2 P04 for details). The main steps of this pro- 
cedure are identification of local density peaks, and analysis 
of the density, circular velocity, and velocity dispersion pro- 
files with simultaneous removal of unbound particles. The 
final profiles, using only bound particles, are used to calcu- 
late the halo properties such as the circular velocity profile 
y circ = [GM(< r)/r] l/2 and the maximum circular velocity, 
Vmax- The completeness limit for halo identification using this 
procedure is ~ 50 particles. The corresponding v max limit is 
~ 80 km s _1 , and the associated mass limit is ~ 5 x 10 10 M o . 
The best-fit power law to the v max - M^q relatiorf^jderived in 
these simulations at z = 1 is 



log M m = 4.57 + 3.21 logv n 



(16) 



where the velocities are in units of km s _1 and masses are in 
units of ]-T l M Q . 

The identified halos were then classified into host halos 
whose centers are not located within any larger virialized sys- 
tems, and subhalos which lie within the virial radius of a 
larger system (see |Tasitsiomi et al.||2004| for details of this 
procedure). Briefly, a halo is classified as a subhalo if its 
center is within r\%§ of the center of a more massive halo, 
where r^o is the radius which corresponds to a mean spheri- 
cal overdensity of 180 relative to the mean density at the given 
redshift. In the real Universe, centers of host halos can be 
identified as locations of the groups' central galaxies, while 
subhalos correspond to satellite galaxies. 

4.2. Model of the AGN population 

Our data can constrain the following two basic properties 
of the AGN population. First, the overall clustering amplitude 
observed in % p (r p ) constrains the mass scale of the AGN host 
dark matter halos. Note that it is the mass scale of host halos, 
not of individual subhalos in which the AGNs might reside, 
which is constrained by the % p {r p ) amplitude. For example, 
a population of small subhalos with v max = 50 km s which 
are located within larger halos with v max = 300 km s _1 has 
the clustering length which closely matches that for larger, 
parent host halos, and not that for the entire population of 
Vmax = 50km s _1 halos. Therefore, the first parameter of the 
model we use to characterize the spatial distribution of AGNs 
is the minimum v max for the halos which can contain the X-ray 
AGNs. The AGN can be located either at the center of such a 
halo or in any of its smaller subhalos. 

Note that we choose to characterize the halos using their 
Vmax rather than the virial mass for the reasons outlined in 



Kravtsov & Kly pin (1999), Nagai & Kravtsov| ( |2005l ), and 
Conroy et al.| ( |200o| ! The maximum circular velocity is a 
more direct measure of the depth of the halo potential wells. 
It reflects the central properties of the halo better, and is less 
subject to the effects of tidal stripping, than the halo virial 
mass which is dominated by the matter at large radii. There- 
fore, we can expect that the stellar content of the halo and 
all baryonic processes in the center, including the AG N activ- 
ity, are better correlated with v max than with M v j r . |Nagai & 



Kravtsov (2005 ) and Conroy et al.| ( |200"6| argue further that 



the best indicator for the stellar mass of subhalos is their 
before accretion onto the host halo. Since our results are not 



14 Hereafter, we use the mass defined within an overdensity threshold of 
180 with respect to the mean density of the Universe at the given redshift. 
The correspondence between v max and mass is quoted for z = 1 unless the 
redshift is stated explicitly. 
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Fig. 8. — The projected correlation functions of halo samples at z = 1, calculated for V m 
equivalent to a uniform scaling of the correlation amplitude at all separations. 




250 and 350km s 1 and / = 0.5. The trend with V ra j n is mainly 



very sensitive to the v max threshold for subhalos (see below), 
we do not make this distinction. 

For a fixed threshold of the host halo v max , the correla- 
tion length of astronomical objects only weakly depends on 
whether these objects are located in the halo central galax- 
ies or in smaller subhalos. As discussed above, the ratio 
^„(n) I p (r p ) should be a much more sensitive and direct in- 
dicator for the satellite fraction. We parameterize this fraction 
by the probability, /, for objects to reside in the central galax- 
ies of host halos; 1 - / is, therefore, the probability for objects 
to reside in any of the smaller subhalos. The two parameters, 
the v max threshold, V m j n , and the probability /, fully specify 
the relation between our model AGN population and the dark 
matter halos and subhalos identified in the numerical simula- 
tions. y m i n is primarily constrained by the correlation length 



of tjpjrp), and / is mostly constrained by the ^„(n)l^ p {r v ) ra- 

ticB 

Algorithmically, we simulate the AGN locations by ran- 
domly drawing the halos and subhalos from the simulation 
box according to the parameters V m i a and /. First, we select 
all bound structures (both halos and subhalos) with v max > 
80km s _I . This threshold is slightly higher than the resolu- 
tion limit in the simulations. We verified that the final results 
are nearly the same when this initial threshold is varied be- 
tween 50km s and 100 km s _1 . We then select only those 



15 Note that we assume that / is independent of v max and we have not 
explored the models in which / is very different for the most massive halos. 
Effectively, our derived constraints correspond to a mean / in the velocity 
range < 

Vmax 5; ^V^n, containing ~ 90% of halos with V m ax > Vmin 
jKlypinetal.|2010^ . 
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subhalos which are contained within halos with v max above the 
given value of V m \ n . We then randomly put a small number of 
objects (10-100) within the simulation bo?|^] With the prob- 
ability 1 - /, the object is associated with one of the selected 
subhalos, and with the probability / it is put in the center of 
one of the halos. This procedure is repeated multiple times 
randomly selecting one of the box axes as the line of sight. 
Using these simulated objects, we derive a model of ^ p (r p ) 
and £r(7r) for each combination of y m ; n and /. 

4.3. Model Correlation Functions i; p (r p ) and £;„(n) 

The correlation functions i; p (r p ) and tj„(n) were derived 
from the simulation outputs on a grid of parameters within 
the range V min 6 [200; 370] km s - ° and / e [0; 1]. Exam- 
ples are shown in Fig.[7jand[8] Comparison of the two panels 
in Fig. 17] illustrates the effect of / on the correlation func- 
tions. For smaller / (more objects located in satellite galaxies, 
left panel), the projected correlation function % p (r p ) shows a 
stronger component at d < 1 hr x Mpc in excess of a power 
law extrapolation from larger separations. This excess is at- 
tributed to the "1-halo" term in the analytic halo model of the 
correlation function. At the same time, t; n {n) shows a larger 
suppression of the correlation amplitude at d < 1 A -1 Mpc 
with respect to £ p (r p ) at the same separations, and a stronger 
enhancement over % p (r p ) at d - 2 — 10 /z~' Mpc for small 
/. This is the consequence of a stronger "finger of God" ef- 
fect in the case when more objects are located in the satellite 
galaxies. Unfortunately, the statistical uncertainties in the real 
data do not allow a detailed modeling of the observed i; p (r p ) 
at small separations. Modeling of the % n (n) at small sepa- 
rations is further complicated b y the effect of uncertainties 
in the redshift measurements (§ |3.2| >. At large separations, 
> 0.1 of the simulation box size, the correlation functions 



derived from the simulations are not reliabl e ( |Kravtsov et al.| 
2004} [Colin et al.|2006||Conroy et al.|2006| . Taking all these 
considerations into account, we will match the observed and 
model correlation functions in the intermediate range of radii, 
d = 1- 12ft" 1 Mpc. 

First, we compute the correlation length, ro, for each com- 
bination (Vmin,/)- This is done by fitting a power law func- 
tion, (r p /ro) _r , to £; p {r p ) in the range r p = 1 — \2lT l Mpc. 
We then compute the ratio ^„{n)/^ p (r p ) (an example is shown 
in Fig. [9]) and fit it in the same range of separations with a 
modified log-normal function, 



4>{d) = 1 + A exp 



(log d - log do) 2 



2D 2 



(17) 



The index g is fixed at the mean best-fit value for all (V m i n , f) 
combinations, g = 1.2. The ^ n (ji)l^ p (r p ) ratio derived from 
the simulations shows substantial variations related to cosmic 
variance (this can be estimated by comparing the results for 
three viewing angles, see errorbars in Fi g. [91 . Therefore, we 
need to smooth the results of the fit by eq.|17| This is achieved 
by fitting low-order polynomials to the parameters D, do, and 
A as a function of V m ; n and /. We found that an adequate 
description is achieved if the best-fit values of D and do are 
approximated as a linear function of /, and A(V m i n , /) is fit by 
a second-order polynomial. An example of the fitting function 
derived from this smooth map is shown in Fig. [9] Due to the 
size of the simulation box, the uncertainties in the smoothed 

16 This is done to approximate the low space density of Chandra/Bootes 
AGNS. 




100 



Fig. 9. — The ratio of the model correlation functions computed for / = 0, 
0.5, 1 , and V m j n = 310 km s _ 1 . Uncertainties at each separation are estimated 
using the variance of the ratio computed for 3 different viewing angles and 
reflect mostly the cosmic variance wit hin the simulation box. The solid lines 
shows the analytic approximation (eq. |17[ . Note that the analytic fits are de- 
rived from the global model with smoothed trends of parameters as functions 
of / and V m j n (see text). The models are shown in the range 1 - 12ft-' Mpc 
which we use in in applying these models to the data. 

model are still finite, but we verified that they are negligible 
compared to those in the data. 

4.4. Application to the Data 

In applying the correlation function model to the data we 
avoid including any sensitivity to the slope of the correlation 
function. The primary motivation is that our method assigning 
the AGN locations to the dark matter halos may be overly sim- 
plistic to correctly predict the details of the correlation func- 
tion shape. Also, the cosmological parameters used in the 
simulation are slightly different from the currently accepted 
values, resulting in a systematic difference in the shape of the 
perturbation power spectrum in the simulated and real uni- 
verses. This said, the models derived from simulations do 
provide a good fit to the % p (r p ) data (see discussion in §5.1 
below). 

Based on these considerations, our^ 2 includes two compo- 
nents. First, we use the value of t he co rrelation length derived 
from fitting the % p (r p ) function (§ 3.2 1. Second, we use the ra- 
tio x = %„{n) 1 1; p (r p ) in the range of separations 1-12 hr x Mpc 
(the data at n < 1 h~ x Mpc are not used because they are likely 



affected by the redshift measurement uncertainties, see 5 3.2 1 
For halos with circular velocities v max ~ 300 km s~' (as indi- 
cated by the amplitude of the AGN correlation function, see 
below), the "fingers of God" extend to v max /H ~ 3 hr l Mpc, 
just in the middle of this range of separations. Formally, the 
constraints on the parameters of the AGN population model, 
y m i n and /, are derived using a^ 2 function computed as 



* 2 = 



(r 



„mod\2 



rO 



z 



(x - x mod ) 2 



(18) 



where the summation in the second term is over the data 
points in the 1 - 10 h l Mpc sepa ratio n range, and the model 
functions are those described in §|4.3| 
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Fig. 10. — Left: Dependence of ro on V m j„ for the objects at the centers of distinct halos (/ = 1) for the simulation outputs at z = 0.5, 1, 2, and 3. There is 
almost no change with redshift for z 55 2.5. Right: the ratio of ^(ff)/^,('"p) for the population model with V m \ n = 310km s _1 and / = 0.5 (50% of objects are in 
the satellite galaxies of the v max > 310 km s _I halos). Any difference in this ratio between different simulation snapshots is within the uncertainties (estimated 
from analyzing 3 different projections for each snapshot). 



4.4.1. Redshift Trends 

The procedure described above provides a correlation func- 
tion model at the redshifts where the simulation outputs were 
saved. The Chandra/Bootes AGNs span a wide redshift range 
and we have to bin the data into several redshift intervals to 
achieve an acceptable level of accuracy for the correlation 
function measurements. Therefore, we need to account for 
any z-dependent trends of the correlation function models. 

Fortunately, for our choice of observables, ro and the 
%„(n) I % p (r v ) ratio, the redshift trends are very weak. This is 
illustrated in Fig. 
V, 



10 



The left panel shows ro as a function of 
for a population model with / = 1 (all objects are at the 
centers of distinct halos) for the simulation outputs at z = 0.5, 
1, 2, and 3. Obviously, there is almost no change in r ( > for a 
fixed V m j n at z < 3. Any changes are much smaller than the 
uncertainties of our ro measurement even for the full sample. 
Therefore, we conclude that the model ro as a function of V m i n 
does not evolve over our redshift range of interestPjThe ratio 
£r(7r)/£ p (fp) also shows little, if any, evolution with redshift. 
The right panel in Fig. 10 shows the results for the population 
model with V min = 310km s _1 and / = 0.5 (50% of objects 
are in the satellite subhalos of v max > 310km s _1 halos). Any 
difference between the simulation outputs is within the uncer- 
tainties (estimated from analyzing three different projections 
for each simulation output). 

The lack of evolution in the clustering model over our red- 
shift interval (and also the lack of detectable evolution of ro 



4.4.2. Adjusting the Results to a Low-o~$ Cosmology 

Finally, we note that the numerical simulations we use were 
run assuming a high value of the power spectrum normaliza- 
tion, cr 8 = 0.9 at z = 0. This results in an incorrect predic- 
tion of the correlation amplitude for halos of a given mass, 
and thus slightly biases the derived parameters of the AGN 
population model, in particular, y m ; n . Obviously, it would 
be best to use the simulations performed for the currently 
favored cosmological model with <x 8 w 0.8 but in general, 
these were unavailable at the time of this investigation. In 
Appendix [A] we describe a procedure which can be used to 
scale the results from the comparison with the simulation to 
any desired cosmology. In particular, if we use the best-fit 
flat ACDM cosmological model derived from the joint anal- 
ysis of the galaxy cluster mass function and other cosmolog- 
ical d atasets, cr 8 = 0.786 and Q. M - 0.268 (Vikhlinin et al. 



with z, see § 5.2 below) indicates that we can safely combine 
the data over the entire redshift range in the sample. Further- 
more, there is no need to weight the models with the redshift 
distribution — one simply can use the results for the simu- 
lation output at z = 1. We take this approach i n fitt ing the 
parameters of the population model, V mm and / (§ |5.1[). In ad- 
dition, we constrain the evolution of V m j n with z (§ |5.2| i under 
the assumption that the fraction of AGNs in the subhalos does 
not evolve (i.e., using a fixed / derived from the analysis of 
the entire sample). 

17 Note that ro as a function of mass does evolve with redshift, as expected. 
However, this evolution appears to be canceled by the evolution in the M - 
v max relation and the trend of ro with M at a given redshift. 



2009 ), the masses reported below should be scaled by a factor 
of 0.69, the V m in's decreased by 10%, and the number density 
of objects increased by 20%. 

5. MODELING RESULTS 

5.1. V m i n and the Satellite Fraction for the Entire Sample 

Figure [TT] shows the combined constraints on the popu- 
lation model parameters, V m ; n and /, obtained from fitting 
the full sample of Chandra/Bootes AGNs. The best-fit ve- 
locity threshold is V min = 320 ± 20 km s" 1 (68% CL one- 
parameter uncertainty). At z = 1, this corresponds to Migo 
4. 1 x 10 12 ft -1 M , or 2.4 x 10 12 h~ l M G after co rrectin g for the 
lower-cr g cosmology (see § 



4.4.2 



and Appendix|Aj | We thus 
conclude, in agreement with the earlier studies (Gilli et al. 



of 



2009 1 [Hickox et al.12009) , that the X-ray AGNs at high red 
shifts reside in small galaxy groups with masses of a factor o 
a few above the present-day mass of the Milky Way. 

As expected from the striking agreement of the observed 
% p (r v ) and %„{n) projected correlation function, the best fit is 



18 [Hickox et al.| {2009] quote a higher mass, ~ 10 B /r' M e for an X- 
ray AUJNs sample with the mean (z) = 0.51, from which they measure 
r = 4.8 /r' Mpc. The main source of the difference, as explained in |Conroy] 



|et al. |(200 8 see their p.l 195), is a more accurate model for the matter power 
spectrum at galactic scales deployed in modern simulations such as those we 
use here. 
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Fig. 1 1 . — Combined constraints on halo occupation parameters, V m ; n and 
/, obtained from the full sample of the Chandra/Bootes AGNs. The best fit 
(Vmia = 320km s and / = 1) is indicated with a star. The contours cor- 
respond to 68% and 95% confidence regions for two interesting parameters 
(AX 2 = 2.3 and 6.2, respectively). The quantity 1 — / represents the fraction 
of AGNs residing in the satellite galaxies. We also indicate the dark halo 
masses corresponding to two values of Vmin- 

/ = 1, i.e. all AGNs are located at the centers of distinct 
dark matter halos. The 90% CL upper limit on the satel- 
lite fraction is 1 — / < 0.12. This is significantly lower 
than the fraction of suitable subhalos in the simulation box. 
For example, within the 1486 halos with v max > 320km s _I 
there are 683 Milky Way-type subhalos (v max = 220km s _1 ), 
and many more smaller galaxies. If all these galaxies had a 
uniform probablity to host an X-ray AGN, we would expect 
l-/>0.32. 

If we consider only the objects with v max > 320 km s _1 , 
11% of them are satellites in more massive halos, which is 
comparable to the upper limit for 1 — / derived from our 
model. However, we still can exclude a uniform probabil- 
ity for hosting an AGN in this case, because the v max > 
320 km s _1 subhalos by definition are members of more mas- 
sive halos, and thus have higher peculiar velocities than most 
of the satellites in our (V m i n , /) model. To illustrate this point, 
we computed the project correlation functions for all halos 



and subhalos with v max > 300 km s (Fig. 12 1. This pop- 
ulation has nearly the same correlation lengthas our best-fit 
model with (Vmin,/) = (320, 1). In the tj p (r p ) function there 
is a strong "1-halo" component at r p < 0.7 Mpc. We do not 
observe such a component in the AGN autocorrelation func- 
tion, even though the quality of £p(r p ) data at small scales 
is clearly insufficient for distinguishing the models. How- 
ever, the v max > 300 km s dark matter subhalos are mainly 
satellites in much more massive halos (median v max for their 
host halo is 530km s _1 ). They have high peculiar velocities 
and show a strong excess in ^(tt) over £ p (r p ) at separations 
1 - 10 /z~' Mpc, not present in the data. 

Therefore, we can conclude that the X-ray AGNs at z — 1 
tend to avoid massive galaxies in the outskirts of yet more 
massive groups and clusters, or satellite galaxies in the v max > 
300km s _1 galaxy groups. Instead, the AGNs are preferen- 
tially located at the centers of distinct dark matter halos. 

5.1.1. Comparison with Previous AGN HOD Studies 

Our conclusion on the preferential location o f AGNs in 
the halo central g alaxies contradicts the results of Padmanab- 
han et al. ( 2009 P09 hereafter) who observe the presence of 



the 1-halo term in the cross-correlation function of optically- 
selected z < 0.6 quasars and Luminous Red Galaxies (LRGs), 
and use this to conclude that a large fraction of the AGNs is 
hosted by satellite galaxies. 

There are obvious differences in the AGN populations stud- 
ied in P09 and our work. A striking mismatch in the derived 
host halo masses — M > 3 x 10 10 h~ x M Q estimated in P09 is 
a factor of ~ 10 lower than our minimum mass — indicates 
that we indeed may be dealing with different object popula- 
tions. However, we also note a significant difference in the 
underlying methods. While our results follow from the AGN 
aMfo-correlation function, the P09 approach is based on crass- 
correlation of the AGNs with LRGs, the tracers typically lo- 
cated in the halos a factor of ~ 30 more massive than the AGN 
hosts (see, e .g., Table Al i n P09). 

Recently, |Miyaji et "aT] f2011| Mil hereafter) published 
the halo occupation analysis for a sample of X-ray selected 
AGNs. Their AGN sample is at lower redshifts (z < 0.36) 
but closely matches ours in terms of the estimated host halo 
mass. Similarly to P09, Mil detect a 1-halo term in the 
cross-correlation function of AGNs and LRGs, qualitatively 
indicating that some fraction of the AGNs may be located in 
the satellite galaxies. However, even if the 1-halo term is de- 
tected, the best-fit AGN HOD models of Ml 1 imply that the 
AGN incidence rate in the satellite galaxies decreases with 
increasing halo mass, approximately in line with our conclu- 
sions above. 

In summary, we can only speculate whether any differences 
in our conclusions with those of P09 and Mil are due to 
not identical AGN samples (optical vs. X-ray selection, or a 
substantial mismatch in the typical X-ray luminosities in the 
Mil and our samples), or different analysis methods (veloc- 
ity space auto-correlation function of AGNs vs. their spatial 
cross-correlation with more massive tracers). A clear answer 
can be provided by detection of, or strong upper limits on, the 
1-halo term in AGN £ p (r v ) awfo-correlation function. Unfor- 
tunately, the present data quality is insufficient for this pur- 
pose. 

5.2. y m ; n as a Function ofRedshift 

Splitting our Chandra/Bootes AGN sample into 4 subsam- 
ples, ~ 320 objects in each, we can measure the correlation 
length, ro, as a function of redshift if we hold the slope of the 
correlation function fixed at the best-fit value obtained for the 
entire sample, y — 1 .84. The ro in the highest redshift bin, 
z = [1.63 - 5], is poorly constrained, while in the other three 
samples we obtain reasonably accurate values, given in Ta- 
ble fl] and shown in Fig. 13 There is almost no change in r ( > 
over the redshift interval z — 0.17 - 1.6, with all the measure- 
ments being consistent with the average ro = 6.41 ft -1 Mpc 
derived for the entire sample. 

Assuming further a fixed / = 1 at each z, as indicated by 
modeling of the redshift-space distortions in the entire sam- 
ple, we can convert the best-fit values of ro at each redshift to 
the threshold circular velocity for the parent dark matter ha- 
los. The results are shown in Fig. 14 We find no detectable 
trend of V m j n or the corresponding mass threshold, Mr, m (see 
Table [TJ, with redshift, either. This appears somewhat coun- 
terintuitive because in a flux-limited sample, such as ours, the 
objects at higher redshift have higher intrinsic luminosities, 
and we might expect them to be located in more massive dark 
matter halos. We note, however, that the studies of optically 
selected QSOs also indicate a wea k or no trend of clustering 
length with the object luminosity ( [Croom et al. 2005 Shen 
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12. — (Left panel:) Model correlation function for all objects in the simu ation output (both halos and subhalos) with v max > 300 km s (dashed line) in 



comparison with the observed f p (r p ). Approximately 11% of objects with such v 
the excess in the correlation power at < 0.7 h 



are satellites in larger halos, which results in the "I -halo" term visible as 
Mpc relative to the model without subhalos (solid line). The peculiar motions of such halos lead to strong 



distortions of £ n (n) with respect to f p (r p ) so that the ratio of <5r0r)/<fn('°p) reaches ; 
(rignt panel) 



2 in the 2 - 6 A 1 Mpc separation range, inconsistent with the observed ratio 



TABLE 1 

Clustering modeling results as a function of redshift 



Redshift range 



£med 

(2) 



(3) 



'•(i 
(4) 



(5) 



(6) 



(7) 



"halo 

(8) 



"/"halo 
(9) 



0.172-0.555 0.374 42.58 5.52 ±0.56 277 ± 28 3.7 ±1.1 1.54 xlO" 4 



i Q9+0.56 



0.116 ±0.035 



0.555 - 1.000 0.738 43.34 7.69 ± 0.99 383 ± 35 8.4 ±2.5 6.12xlQ- j 



50 +G:2i 
V JV "II. 



0.123 ± 0.037 



1.000- 1.630. 



T779" 43~S5 6.89± 1.11 355 ± 64 4.9 ± 2.7 2.56 x 10~ J 



"■^ y -0.17 



0.052 ± 0.029 



Note. — Column (3) — mean X-ray luminosity in the 0.5-2 keV band (rest-frame). Column (4) — best-fit correlation 
length assuming a fixed slope of y = 1.97, in units of h comoving Mpc. Column (5) — threshold maximum circular velocity 
of the host dark matter halos, in units of km s - ' (j |5.2| . Column (6) — the halo virial mass corresponding to V m j n , in units of 
10 12 /i 1 M Q . Column (7) — comoving number density of X-ray sources at z me d, in units of /i 3 Mpc 3 . Column (8) — comoving 
number density of dark matter halos with v max > V m i n at Zmed, in units of 10~ 3 /i 3 Mpc -3 . Column (9) — probability for the 
Chandra/Bootes AGNs to be in active state computed as the ratio of number density of AGNs and the dark matter halos with 
max > ^min- 
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Fig. 13. — The correlation lengths, ro, measured in 3 redshift intervals, z = 
[0.172 - 0.555], [0.555 - 1.000], and [1.000 - 1.630], containing 320, 307, 
and 344 objects, respectively. The values of ro were obtained by fitting the 
projected correlation function, w p (r p ), in each bin assuming the fixed slope 
of the correlation function y = 1.84 (the best-fit value for the entire sample). 



et al. 
etaT 



2009). Clustering analysis of the SDSS quasars (Ross 
2009 ) shows mild or no evolution of the real-space cor- 



relation length at z < 2. 



Fig. 14. — Minimum rotational velo city for dark matter halos, corresponding 
to the correlation lengths in Fig. |13| 

5.3. Constraints on AGN Duty Cycle 

Having estimated the mass scale and therefore the space 
density of Chandra/Bootes AGNs, we can formally compute 
their duty cycle following the approach of Martini & Wein- 
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berg (2001 1. Under the simplifying assumption that the halo 
lifetime is approximately independent of mass, the probability 
for the AGN to be active is simply 



«AGN 
"halo 



(19) 



where «agn is obtained from the dN/dz fit (eq.j2j for the me- 
dian redshift for the given subsample. The limiting X-ray lu- 
minosity is ill-defined for our sample because the X-ray de- 
tections extend to very low limits in terms of the number of 
detected X-ray photons (> 4, see Rent er et al.|2005| , and in 
this regime, a wide range of possible intrinsic intensities cor- 
responds to the given number of photons (Renter & Murray 
2003). However, as a guide for the typical luminosity scale 
one can use the median log L x reported in Table [T] With these 
caveats, the probabilities given by eq. 19 and reported in col- 
umn (9) of Table [T]correspond to the probability for the dark 
matter halos with v max > V^n to host an AGN with the in- 
stantaneous soft-band X-ray luminosity of order (log L x ) or 
above. These P act i V e generally decline with mean luminos- 
ity and/or redshift. Interestingly, our results indicate that 
AGNs are quite common at low redshifts — approximately 
10% of dark matter halos with v max > 280 km s (or mass 
Migo > 3.7 x 10 12 hr l M ) host an AGN with a soft-band X- 
ray luminosity of ~ 4 x 10 42 erg s _1 or above. 

6. DISCUSSION AND CONCLUSIONS 

We measured the clustering properties of X-ray selected 
AGNs using a sample of 1282 sources with spectroscopic red- 
shifts from the 9.3 deg 2 Chandra survey in the Bootes region 
— one of the most accurate such measurement to date. In 
agreement with previous studies of X-ray and optically se- 
lected AGNs, we find that the real-space correlation function 
can be approximated in the radial range 1 - 20/i _1 comoving 
Mpc by a power law with a slope of y — 1 .84 + 0.12. The cor- 
relation length is w 6.4 h Mpc, showing only weak trends 
with redshift at z — 0.17 -1.6 (or with X-ray luminosity). 

Matching the observed clustering properties of the Chan- 
dra/Bootes AGNs to those of dark matter halos in the high- 
resolution cosmological simulations, we find that the X-ray 
AGNs reside in halos with the maximum rotational velocity 
» 320 km s _1 , or with total masses ~ 4.1 x 10 12 /z -1 M o , also 
with no detectable redshift trend. The lack of a redshift or 
luminosity dependence of the AGN clustering is inconsistent 
with the scenarios in which the AGN luminosities in the active 
state are similar fractions of the Eddington luminosity. How- 



ever, it can be explained in the scenario in which the AGN 
activity is triggered by major mergers of gas-rich galaxies, 
and the instantaneou s luminosi ty passes through many levels 
after each trigger (Hopki ns et al.|2006| . 

Our results reveal another interesting aspect of the AGN 
clustering which was predicted in Hopkins et al. ( 2008[ ). The 
redshift measurements in our sample are sufficiently accu- 
rate to detect peculiar motions of objects in excess of ~ 
100km s _1 . The comparison of the two-point correlation 
functions projected on the line of sight and on the sky plane 
reveals no signatures of the redshift-space distortions, which 
allows us to put limits on the fraction of AGNs located in the 
satellite subhalos within the host dark matter halos. We find 
that the X-ray AGNs are predominantly located in the central 
galaxies of the host dark matter halos and tend to avoid satel- 
lite galaxies. Quantitatively, we limit the fraction of AGNs 
in non-central galaxies to be < 0.12 at the 90% CL. We also 
exclude the model in which the probability for a galaxy with 
Vmax > 300km s -1 to host an AGN is the same for central 



galaxies and satellite galaxies in more massive halos (§5.1 



The central locations of the quasar host galaxies are expected 
in the trigger model because mergers of equally-sized galaxies 
preferentially oc cur at the centers of dark matter halos ( |Hop-| 
|kins et al.|2008) . 

Finally, we compared the number densities of the Chan- 
Jra/Bootes AGNs to that of the dark matter halos with the 
mass corresponding to the AGN clustering amplitude. We find 
that the fraction of halos with active X-ray AGNs decreases 
with increasing z — and, correspondingly, with L x — proba- 
bly reflecting a lower probability for an object to have a higher 
instantaneous luminosity. At the lowest redshifts in our sam- 
ple, Chandra probes such low luminosity that X-ray AGNs 
become quite common. At z = 0.37, the C/ianafra-detected 
sources are located in more than 10% of the dark matter halos 
with v max > 280km s" 1 or M > 3.7 x 10 12 /T 1 M„. 
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APPENDIX 

COSMOLOGY-DEPENDENT SCALINGS OF THE MODEL CORRELATION FUNCTION 

Using the simulation outputs directly has many advantages over the analytic HOD model in testing different models of pop- 
ulating the dark matter halos with astronomical sources. However, there is an important disadvantage. The simulations are 
performed for a certain combination of the cosmological parameters, which may or may not match the currently favored cosmo- 
logical model. For example, the numerical simulations we use in this work were performed assuming a high value of the power 
spectrum normalization, cr 8 = 0.9 at z = 0. This results in a slightly incorrect prediction of the correlation amplitude for halos of 
a given mass, and thus slightly biases the derived parameters of the AGN population model. 

Fortunately, as we discuss below, the simulation-based models can be easily adjusted to the "correct" cosmology by simply 



rescaling the halo masses by a small factor. A more elaborate procedure is described by Angulo & White ( 2010[ ) who show that 
even the raw simulation outputs can be rescaled to the correct cosmology. 

Basics 

The correlation function of objects with mass M at a given redshift, z, can be written as 

{(r,M) = b 2 (M){ DM (r), (Al) 

where b(M) is the bias factor for halos with mass M, and £dm(/) is the correlation function of all dark matter particles. Since 
we consider the correlation function measurements at the separations bracketing 8/i~' Mpc, we can assume that with sufficient 
accuracy, ^dmM k o~l(z) = ^loDiz) 2 , where D(z) , is the linear perturbations growth factor for the given cosmology. 

The bias factor, b(M), is, in turn, a function of the linear perturbation amplitude at the mass scale M, b(M) — b(cr(M)) (Kaiser 
1984 Mo & White|1996 Sheth & Tormen|19 99 ). Therefore, we conclude that the dependence of the correlation function model 



on cosmology is through erg, D(z), and cr(M). 



Scaling for o~{M) 



The scales corresponding to the galaxy-sized objects are sufficiently small that the slope of the ma tter power spectrum and 
hence the mass dependence in cr(M) is insensitive to the underlying cosmology (Bard een et al.|1986 i. In a fairly broad range 



of parameters around t he "concordance" cosmological model, we find that cr(M) computed with the full transfer function (e.g., 



Eisenstein & Hu|l998| can be approximated as cr(M) — Ax M 016 in the mass range 10 12 - 10 13 M Q . The amplitude, A, of this 
power law approximation scales with cosmology as 

A oc cr 8 x D(z) x s(£l M , h, ...). (A2) 

In this equation, the factor s{Q.m, h, ...) represents the correction for the diff erent power spectrum shapes in different cosmologies 
(it is mostly a function of the product fi M /!, see Bond & Efstathiou 1984). We neglected such a factor in %dm(/) because we 
consider the cosmological models in which the perturbations are normalized by erg, the fluctuation amplitude at approximately the 
midpoint of the observed range. However, the scale corresponding to galaxy-sized objects for which cr(M) has to be computed, is 
sufficiently far from 8 hr l Mpc so that the correction can be important. For example, we find that s\ j s^ — 1 .03 for the flat ACDM 
cosmologies with {Q-u, h) = (0.3, 0.7) and (0.268, 0.71) (more power in the high-0 M case). 
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Fig. 15. — Cosmological functions <x(M) and b(M) computed for two different fiat ACDM cosmologies, <xg) = (0.3, 0.9) and (0.268, 0.786). The left panel 
shows cr(M) computed at z = 1 (solid lines) and the power law approximations, <x oc M~° 16 (dotted lines). In the right panel, the linear bias, b(M), is plotted as a 
function of cr(M) for the same two cosmologies and for z = 1 and 2 (the mass range is M= 10 12 - 10 13 as in the left panel). 



Scaling for b(cr(M)) 

The key notion for our computation is that the halos with M - 10 12 - 10 13 are "rare" objects at z > 1 (cr(M) 



high bias factors strongly dependent on cr(M). The bias as a function of cr(M) can be computed using the Sheth & Tormen 
approximation: 

aS 2 c /cr 2 -l 2p 



1 or less) with 
T999) 



1 + 



6 c [l+(a6 2 c /o- 2 f] 



(A3) 



with parameters a = 0.75 and p = 0.3 (adopted from Hu & Kravtsov|2003| see their page 704), and 6 C = 1 .69 is the threshold 

nated universe. The bias computed from 



for spherical collapse in a matter-dominated 
different redshifts is shown in the right panel of Fig 

b oc <x~ L6 at z = 2 and b oc cr~ 11 at z = 1. 
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computed 
For the mass range 



this equation for the two cosmologies and two 
10 12 - 10 13 M , the bias can be approximated as 



in our mass range. 
-1.6...-U For 



Scaling of Mass Derived from the Two-Point Correlation Functions 
Since <x(M) is a weak function of mass, the linear bias is also a weak function of M, b(M) oc m°- 2 -0- 2 5 
However, cr(M) scales linearly with cr g , and the bias also shows a strong dependence on this parameter, b oc cr g 

the correlation function of halos, we have £ = £dm^ 2 00 °"g 12 ■ ■ -^s 04 ' therefore, arrive at a somewhat counterintuitive 
conclusion — the amplitude of the correlation function of galaxy-sized objects at high z is lower for models with high values of 
cr 8 because of a strong dependence of the linear bias on the underlying amplitude of the density perturbations, cr(M). 
When we determine the mass scale of objects from their correlation function, we effectively solve the equation 



£ DM b\M) oc o-\ D(zY b\M) = C, 
where C is a constant provided by the data. Inserting the scalings derived above, we can rewrite this as 

o-\D(zf [<t$D(z)sM-° m ]~ 



-2a 



c 



(A4) 



(A5) 



where we assume that b(M) is approximated as a power law of cr(M), b oc cr a (e.g., a 
From here we have, 

M 0M oc[o- & D(z)f- 1 s a 



-1.6 for our mass range and z — 2). 

(A6) 



If we scale the object masses by this equation, the halo clustering properties in the two different cosmologies will be very similar. 

As an example, consider the scaling betwee n the cosmology used in the simulations, Q.m = 0.3, h = 0.7, <r g = 0.9, and the 
galaxy cluster-normalized best-fit model from | Vikhhnin et al. (2009), Q M = 0.268, h = 0.715, cr g = 0.786. At z = 1.03, the 
median redshift of our Chandra/Bootes X-ray sample of AGNs, D(z) = 0.616 and 0.604 for Dm = 0.268 and 0.3, respectively. 
The shape factors, s, are 1 and 1.03 for Q. M = 0.268 and 0.3, respectively. Inserting into eq. A6 we find that the masses which 



would be derived in the D.m - 0.268, cr 8 = 0.786 cosmology are a factor of 1.45 lower than those derived using the correlation 
function models obtained directly from the simulation outputs. 

At z — 2, D(z) = 0.284 and 0.277 for Q.m = 0.268 and 0.3, respectively. Using eq. | |A6| , we find that the object masses should 
be scaled down by factors 1.7 ... 2 for a — 1.4 ... 1.6 (the range of slopes for the z - 2 objects, see Fig. 15 i. 
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Scaling for the Number Density of Objects 

A related question is how to scale the number density of objects with the given clustering properties. In many situations, it 
should be possible to use a mass function model (e.g., the one from Sheth & Tormen||1999| l for such a scaling. However, for 
g alaxy-sized ob j ects i t is possible to derive a simple scaling for the number density. 

|Jenkins et al. ( |2001| l show that the mass function of halos can be written in a standard form, 



dn 1 diner 1 

dM M- d\nM Y 



(A7) 



where B is a cons tant, cr(M) is the linear perturbations amplitude on scale M, and (f(cr) is a "universal" function. For the Sheth 
|& Tormen| ( fT999l l model, 

6 r (ad 2 



2a 

<p(o-) = 032 J — 
V n 



1 + 



— exp -r— ~ 

cr \ 2 cr 1 



(A8) 



where the parameters a and p are the same as in the expression for the linear bias (eq. (A3)). For cr < 1, as is the case for our 
galaxy-sized objects, the Sheth & Tormen / is a slowly varying function of cr which can be approximated as 

<p(cr) * 0.65 +lncr (A9) 

Inserting this into eq. |A7| and taking into account that cr(M) oc M~°- 16 (see above), we have for the number density of objects 
above a mass threshold Mi; m , 



N 



r°° o.i6 
= fix — - 



[0.65 + In cr lim - 0.16 ln(M/M lim )] dM <x MZl (0.49 + In cr llm ) 



(A10) 



where cr lim = cr(M\i m ). Now we can take into account that the limiting mass is found from the conditi on|A4| (so that the clustering 
properties match those observed). If bias is approximated as a power law of cr, b(M) oc cr{MY a (Fig. [If 

<x[cr 8 D(z)] 1/ff . 



(All) 



Equations |A10| |ATT| and |A6| provide the scaling for the number density of isolated halos whose clustering properties match 
the observations. Assuming that the relative number density of satellite and main halos is not very sensitive to the underlying 
cosmology, the same scaling can be applied to any other type of halos. 



In summary, we suggest the following procedure. One takes the simulation outputs and adjusts the mass threshold, My lm so that 
the correlation function amplitude derived for the halos with mass above this threshold matches that observed. One then meas ures 
the number density, N, of su ch halos in the s imu lation box. To convert these quantities to the desired cosmology, one uses eq. |A6 
to scale My lm and equations A10 All and A6 for the number density. For example, we found for the Chandra /Bootes AGNs 
^min = 320km s~', which corresponds to M\ lm = 4.1 x 10 12 /z -1 M . This corresponds to cr lim = 1.25 at z — 1. To adjust these 
results from the (QmiCs) = (0.3,0.9) to (0.268,0.786) cosmology, we need to scale M\ lm by a factor of 0.69, and increase the 



number density by a factor of 1 .2 (the M^l factor is partially compensated for by the change in In cr^ 
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